Rigorous Derivation of the Gross-Pitaevskii Equation 
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The time dependent Gross-Pitaevskii equation describes the dynamics of initially trapped 
Bose-Einstein condensates. We present a rigorous proof of this fact starting from a many-body 
bosonic Schrodinger equation with a short scale repulsive interaction in the dilute limit. Our proof 
shows the persistence of an explicit short scale correlation structure in the condensate. 
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Recent experiments on trapped Bose gases in the con- 
densate phase P, [3] have revived interest in the rigorous 
justification of the Gross-Pitaevskii (GP) theory from 
interacting many-body Hamiltonians. In these experi- 
ments, the gas is initially trapped by a strong magnetic 
field and cooled down at very low temperatures, so that 
almost all particles condensate in the same collective 
mode described by a one-particle wave function, <p(r), 
r e M 3 . Then the traps are instantaneously removed 
and the evolution of the condensate is observed. To de- 
scribe this evolution, the GP theory 0, 0, [B| postulates 
that many-body effects can be compressed into a non- 
linear on-site interaction. In units where H = 1 and the 
mass of the bosons m = 1/2, the GP equation is given 

by 



id t (p t = -A<p t + o-\tp t \ 2 <p t 



(1) 



with the normalization J" |c^(r)| 2 d 3 r = 1. The coupling 
constant is a = 8?rNa where N is the number of particles 
and a is the scattering length of the interaction potential. 

A single orbital theory might indicate that the full N- 
body wave function essentially factorizes. This, however, 
is usually not the case: very short scale interactions may 
introduce strong pair-correlations that substantially in- 
fluence the energy of the system and the dynamics even 
on larger scales. In particular, the emergence of the scat- 
tering length in a is a correlation effect and it is remark- 
able that this correlation can be consistently modeled by 
the same coupling constant along the whole evolution. 

The GP equation thus implicitly postulates a persis- 
tent two-scale structure of the state. On short scales, the 
state exhibits a time independent pair-correlation, while 
on large scales it is given by the product of N copies of a 
time dependent orbital. It is a fundamental problem in 
many-body theory to show that the GP postulation of a 
two-scale structure is a rigorous consequence of the first 
principle many-body Schrodinger equation. While many 
theoretical work addressed this question for low energy 
states, there is no rigorous result from a dynamical point 
of view. In this letter, we shall show that not only the 
two-scale structure is preserved, but it even emerges dy- 
namically for a class of initial data. 

On the mathematically rigorous level, the GP the- 
ory has only been verified for the ground state of the 
trapped Bose gas. Lieb, Seiringer and Yngvason [|| con- 
sidered a model where the scattering length a of the in- 
teraction varies with N, so that Na remains constant. 



For repulsive interactions they proved that the ground 
state energy per particle in the GP limit (N — > oo, 
ao = Na = fixed) is given by a variational principle, 
min{£ G pM : IMI = 1}, where 



(2) 



£ GP (tp) = J d 3 r(|V^(r)| 2 + y cxt (r)|(p(r)| 5 
+47ra |^(r)| 4 ) 



is the GP energy functional. This result was extended to 
rotating Bose systems that model superfluidity in 0] • In 
[H, Lieb and Seiringer have also rigorously established 
that the ground state of the trapped Bose gas exhibits a 
complete Bose-Einstein condensation (BEC) by showing 
that the one particle density matrix of the ground state 
converges to the one dimensional projection |</?gp)(v 5 Gp|, 
where ifcp is the minimizer of the functional £gp- 

The ground state has been known to have a non-trivial 
pair correlation. Mathematically, this was first exploited 
for the homogeneous dilute Bose gas, when Dyson § 
used a specific nearest neighbor correlation in his trial 
state to verify the upper bound on the ground state en- 
ergy obtained from Bogoliubov theory [10( . The proof of 
the lower bound, achieved by Lieb and Yngvason forty 
years later , demonstrates that low energy states nec- 
essarily contain a characteristic short scale pair correla- 
tion. In the GP limit, as shown in [fj, the trapping 
potential introduces a second length scale but it leaves 
the short scale structure intact. 

The works mentioned above justify the two-scale hy- 
pothesis of the GP theory based on the energy minimiz- 
ing principle. In this letter, we outline a mathemati- 
cal proof (see [HI, for details) of the fact that the 
GP theory also describes the dynamics of trapped Bose- 
Einstein condensates. Assuming that at time t — we 
have a condensate, we show, under some conditions on 
the interaction potential and on the initial state, that 
the system still exhibits complete condensation at times 
t 0, and that the condensate wave function evolves 
according to the GP equation {]]). Our results cover 
two types of initial data: type I states with the phys- 
ical two-scale structure as in jf| Q and type II states 
with no short scale structure. For type I initial data, 
the main difficulty is to demonstrate the persistence of 
the two-scale structure even for states that are far from 
the ground state. Energy conservation alone does not 
prevent the destruction of the short scale structure. The 
key observation is that higher moments of the energy are 
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also conserved and they impose stringent restrictions on 
the state. In Theorem 2 below we show that a finite 
second moment of the energy per particle forces a pair 
correlation described by the zero energy scattering func- 
tion of the interaction potential. For type II initial data, 
our results imply that the short scale structure emerges 
dynamically. See the discussion after Theorem 1. 

We now define the model precisely. The Hamiltonian 
for N bosons in three dimensions is given by 



H 



N 



N 
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-Ar 3 + Kxt^-)) + V "( r * - r,-) . (3) 



The trap potential V^ x t(r) is a positive function with 
Ksxt(r) oo as |r — > oo. The interaction potential 
is scaled as Vjv(r) = N 2 V(Nr), where V is a positive, 
spherically symmetric, compactly supported, smooth po- 
tential with scattering length oq. By scaling, the scatter- 
ing length of Vn is a /N. We assume that the strength 
of the interaction, measured by the dimensionless quan- 
tity a = J F(r)|r| _1 dr + sup reR3 r 2 V(r), is sufficiently 
small. 

Let ijjjy be the ground state of (J3]) and let yffi be its k- 
particle marginal density (1 < k < N) with a consistent 



normalization Tr 7^ = 1 . It was proven in [8j that 
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for any fixed k. In particular, ((4]) for k = 1 means 
that ipN exhibits complete condensation. After instan- 
taneously removing the trap, the evolution of the system 
is generated by the Hamiltonian 
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-A r . + Y J V N {r i -r j ). 



(5) 



Theorem 1. Let ipN,t be the solution to the 
Schrodinger equation idtipN,t = HifipN,t with the ini- 
tial condition V>iV,o — ipN) an d ^"/n\ be its one particle 
marginal density. Then, for any fixed time t 6 K, ipN,t 
exhibits complete Bose-Einstein condensation, that is 



IN,t 



\<Pt){yt 



as N 



00, 



where tp t solves the Gross- Pitaevskii equation 
id t <pt = -A(£> t + 8na \<pt\ 2 <Pt 



(6) 



(7) 



with initial data tpt=o = <Pgp- The convergence in (OJ) 
is in the sense that Tr K (7^ — |yt)(¥>t|) for any 
compact operator K on L 2 (R 3 ). 

Our proof, in fact, allows for more general initial data. 
Let H 1 (R 3 ) := {<p(r) : / \<p\ 2 + \V<p\ 2 < 00}. The tech- 
nical condition we need to conclude ©, (J7J is that the 
initial state ipN.o asymptotically factorizes in the sense 
that there exists <p € 7? 1 (R 3 ) and for every k > 1 there 



exists a family k ^}N>k 



with ^- k) e L 2 



and \\(p\\ = = 1 such that 



lim 
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N,0 



= 



(8) 



The function (p is then the initial data for ([7]). As we 
prove in Appendix C of [ijj], the ground state i^n of 
the Hamiltonian ([3]) satisfies J5]). This shows that the 
pair correlations of ip n , which affect the dynamics in a 
non-trivial way, are not in conflict with the asymptotic 
factorization as long as the convergence in (JSJ) is in a 
weaker topology than the energy norm. The local corre- 
lation structure of ipN lives on a length scale 1/N and its 
effect, measured in i 2 -sense in ©, is negligible. How- 
ever, the energy of tj)N does not converge to the energy 
of the factorized state, and, similarly, the convergence in 
((H) does not hold in the H 1 space or in energy sense. 

A typical type II initial data is given by a product 
wave function tpNfi = <P® N with ip 6 H 1 (M. 3 ). For 4>n,o 
the condition ((5|) holds trivially. Our result thus shows 
that a product initial state also evolves according to the 
GP equation ([7]) (with initial data tpt=o = f) despite 
that its total energy is larger than £gp(v)- In fact, it 
is easy to see that the energy per particle in the factor- 
ized state, (tp® N , (Hn /N)ip® N ), is asymptotically given 
by a functional similar to ^ but with 87rao replaced 
with its Born approximation, 60 = / V(r)d 3 r. Note that 
60 > 87rao. This indicates that the excess energy orig- 
inating from the emergence of the scattering length is 
dispersed into modes that do not influence the evolution 
of the condensate. Hence, for product initial data, the 
scattering length emerges dynamically and the energy 
functional of the single particle orbitals does not predict 
the evolution equation. 

Outline of the proof. For k — 1, . . . , N, let 7^ denote 
the fe-particle marginal density of ipN,t with a consistent 
normalization Tr Jn\ = 1. The time evolution is gov- 
erned by a hierarchy of N coupled equations, commonly 
known as the BBGKY hierarchy 



0=1 i<j 
k 

+ (A - fc) £ Tr fc+1 [V N (rj - r k+1 ), 7 ^ 1) ] 
0=1 



(9) 



where Ttu+i denotes the partial trace over r^+i. By the 
Banach-Alaouglu theorem, 7^ has at least one weak* 

limit point, 7^\, in the space of trace class operators. 
Since formally NV N (r) = N 3 V(Nr) -» b S(r), one may 
naively expect from © the time evolution of 7^ to be 
described by the infinite hierarchy 



0=1 



(10) 



+ a V Tr fc+ i [5(Tj ~ r fc+1 ), 



(k+i) 



0=1 



with a = &o- However, the coupling constant in the 

(k) 

correct equations for 7^,^ is a — Sirao- The reason is 

that 7^^ develops a short scale structure which lives 
on the same length scale £ ~ 1/N as the potential Vn. 
Since Jn^ 1 ^ is not constant on this length scale, one 
cannot apply the formal limit NVn — > boS in ^j. 
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In the following Theorem 2 we show that for Tj close 
to rfe + i, the short scale structure of 7^1* is described 
by the function /jvfo — r fc+i)j where /at is the solution 
to the zero energy scattering equation 



-A + ^(r)) Mr! 







(11) 



with boundary condition /zv( r ) ~ * 1 as |r| — > oo. 
By scaling, /^(r) = f(Nr), where / is a solution of 
(— A + ~V) f = with the same boundary condition. 
Therefore the correct value of a in (fTU|) is given by 



d 3 vNV N (v)f N (v) = / d 3 r V(r)/(r) = 8™ ■ 



Note that the short scale structure disappears from 7^ t 
after taking the weak limit, but it still affects the limiting 
macroscopic dynamics. 

The infinite hierarchy (fTU|) has a factorized solution: 
the family jf® t = \(p t ) (tp t \® k , k > 1, is a solution to (fTT)|) 
with a — 87rao if and only if (p t is a solution to the GP 
equation ((TJ. Therefore, to identify the limiting density, 
it suffices to show that (i) any limit point {7^\}fc>i of 

the family of densities {"f]!j\}k=i 1S a solution of the infi- 
nite hierarchy (|I0[) with a — 87rao, and (ii) the solution 
to the infinite hierarchy is unique. The proofs of (i), (ii) 
rely on estimates on the energy distribution of the initial 
wave function and on the factorization property ([5]). 

(i) Convergence to the infinite hierarchy. The follow- 
ing key theorem identifies the short scale structure as- 
suming a bound on the second moment of the energy. 
This is the main ingredient in the proof of (jTUJ) 

Theorem 2. Suppose a is small enough. Then there 
exists a constant C > such that 



(iP, H%i,) > CN 2 



dR 



V rj V, 



^(R) 



(12) 



1 f N {Vi-Yj) 

for all 1 < i < j < N and all tp E L 2 (R 3JV ,dR) 
symmetric with respect to permutations. Here R = 
(ri,.,r,)6l M 

From this theorem, we obtain the a-priori bound 

2 



dR 



V^(R) 



<C, Vi^j (13) 



j /iv( r i - r 

uniformly in N and t provided the initial wave function 
ipN,o at time t = satisfies (4>n,o, -ffjvV'iv.o) < CN 2 . 
Note that / di^dr,- 1 V ri V rj f N {r t - r ^ ) | 2 ~ N thus (H2J) 
identifies the short scale structure of ipN.t with a preci- 
sion 1/iV. Furthermore, since the left side of (12) is a 
constant of motion, (fT3|) shows that the separation be- 
tween the singular short scale structure and the regular 
part of ?pN,t is preserved by the time evolution. 

Outline of the proof of Theorem 2: we define h m = 

~ A r m + \ En^rn V N^n - T m ) for TO = I, ... , N. Then 

H]s = X)m=i and, using the permutation symmetry 
of ip and V > 0, we obtain, for arbitrary i =7= j, 



N(N - 1) 



> 



> {ip,hihjip) 

V N (r t - Tj) 



)(-A rj + 



V N {Yj 



-)i> 



From the definition (fTTj) of /jv(r), we find 

1 v 

(- A ri + ~V N (Ti - Tj))lj} = f N {Ti - Y^U—— 

z JN \ r i — r j ) 

with Li = — A r; +2V ri (log/Ar(ri — r 3 -))V Pj and an anal- 
ogous identity for Lj. Using integration by parts, and 
with </>ij(R) = ip(R)/f N (ri - r 3 -), 

> JdnfUn-Y^v^^m 2 

+ J dR (/|rV 2 log/^) (rj-r^V^.fRJV^^CR). 

Eq. (fT2|) now follows because |V 2 log/Af(ri — r 3 )| < 
Ca |r^ — r 1 — 2 and therefore the second term on the 
right hand side of the last equation can be controlled 
by the first one (for a small enough) using the operator 
inequality r|~ 2 < — CA r (Hardy inequality). 

(ii) Uniqueness of the infinite hierarchy. The first step 
is to prove a-priori bounds in a certain Sobolev norm. 

Theorem 3. Let be any weak limit point of 
then the following estimate holds uniformly in time 



| 7 S t |* := Tr (1 - A ri ) . . . (1 - A r J 7 £ t < C. (14) 
Idea of the proof. By conservation of iT^r along the evo- 



(fe) 



lution, Tvm 



A control on TrH k ~ {k) 



l NT N ,t = LtU n1n,o- a control on irii^^ 
can be obtained through ([8]). The difficulty is that the 
norm ||| • |||fe cannot be directly controlled by TrH^(-) and 
in fact | 7 jy jlfe — ► 00 as N — ► 00 because of the singular 
short scale structure. It is only after taking the weak 
limit N — ► 00 that the short scale structure disappears 
and (fT"4"| can be proven. 

For illustration, consider the case k = 2 discussed in 
Theorem 2. The estimate (fT3|) implies that ^>jv,j(R) ~ 
/jv(rj — r 3 )0t(R) where ^ is a smooth function in 
the variable r; — r 3 . Together with the fact that 
J drjdrj|V ri V Fj . /jv(i"i — ij)| 2 ~ iV, we have that 

JdR j V ri V rj f/'Af(R)| 2 — > 00 as N — ► 00. However, 
since /at — > 1 weakly in L 2 -sense (but not in energy 
sense), a Sobolev estimate on the limit of the density 
matrix of ^j\r,t can be deduced from (jT3|) . 

The uniqueness holds in the Sobolev norm ||| • |||/.: 
Theorem 4. Given a family of densities T = 
{7^}fc>i such that W^f^Wk < C k , there exists at most 
one solution T t — {7t' C ''}fc>i to with r t= o = T and 
such that | 7t \\k < C fe holds uniformly in t. 

Outline of the proof. Iterating the integral form of 



T0[) . we obtain a Dyson series 



It 



(k) 



n-l 



where W t (fe) 7 « = e l * ^J=i ^ 7 (fc) e -« E, fc =1 A, ( and 



(fc) 



dsi . 
Jo 



(16) 



x ^(fe+m-l)^(fe+m)^,(fe+m)^ 
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(k) 
m,t 
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Graphs 



(k+m) 



I 



( k) 

1: Diagrammatic expansion for u) m t as sum of Feynman 



FIG 

graphs with m vertices, 2{k 



leaves and 2k roots 



with the collision operator given by 



-ia Tr 

„(*) 



fe+i 



rfc+i),7 



(M-i)i 



,(fc) 



„ 4 has the same form as uj n t 

(fc+n) 



with 



The error term r/„ 
fc+«) 7 (fc+n) repla( 

To prove the convergence of the expansion (fT5|) as 



wi^ + ™'7^ fc+ "^ replaced by the full evolution 7, 



n — ► 00, we expand each term into a sum of contribu- 
tions associated with certain Feynman graphs. A typical 
graph A contributing to u>^\ is drawn in Fig. [TJ It has m 
four-valent vertices and 2fc + 3m lines. The external lines 
on the left (called roots) correspond to the 2k momenta 
variables of the operator kernel of The 2{k + m) 

external lines on the right (called leaves) represent the 
kernel of 7 ( fe+m ) . The graphical structure of A encodes 
the collision history given in (| 1 6ft . Every line e of A car- 
ries a regularized free propagator, (a e — Pg ± i/i) -1 with 
a momentum variable p e 6 R 3 and a frequency a e G M. 
At each vertex, there is a p- and a a-delta function due 

to momentum and energy conservation. The kernel of 

(k) 

oj m t (whose variables correspond to the 2k momenta of 
the roots) is computed by performing 3m momentum 
integrals and 2k + 3m frequency integrals in each A. 

Because of the singularity of the interaction at r = 0, 
each graph is potentially ultraviolet divergent. Power 
counting suggests, however, that the integrals are fi- 
nite. Suppose we cutoff all momentum integrals at 
p ~ (3 ^> 1, and all frequency integrals at \a\ ~ /3 2 ; 
then the integration volume scales as ( (3 3 ( 3m )+ 2 ( 2 ' c + 3m ) = 
/3 4fc+15m . The m momentum and frequency delta func- 
tions scale as (3~ 5m and the 2k + 3m propagators scale 
as j g-2(2fc+3m) _ a _p r j or j es ti ma te (fT4|) shows that 

7 (fc+m) < /3 -5(fc+m)_ Since 4fc+15m < 5m+2(2fc+3m) + 
5(fc + m), the integrals should be convergent in the ultra- 
violet regime /3 > 1. To make this argument rigorous, 
we use an integration scheme, dictated by the structure 
of the graph. We start by integrating the momenta and 
frequency of the leaves; these integrals are convergent be- 
cause of the a-priori estimates (fl4"|) and they also provide 



a momentum decay on the lines adjacent to the leaves. 
We then iterate this procedure, integrating all momenta 
and frequencies by moving from the right to the left of 
the graph and transferring a suitable momentum decay. 

Conclusion. We have proven that Bose-Einstein con- 
densates evolve according to the Gross-Pitaevskii equa- 
tion. This provides a mathematical description of recent 
experiments on the evolution of initially trapped conden- 
sates. On the theoretical level, our result for factorized 
initial wave functions is more surprising. The emergence 
of the scattering length in the GP equation is a conse- 
quence of the short scale correlations. We show, how- 
ever, that the GP equation is correct even if the initial 
state is uncorrelated (product). Our result thus suggests 
that the many body dynamics builds up correlations on 
the length scale I ~ 1/N in a very short time. Since 
the emergence of correlations reduces the local energy, 
one may ask what happens with the excess energy. Al- 
though we are not able to give a mathematically rigorous 
answer, we believe that the excess energy is transferred 
to incoherent excitations living on intermediate length 
scales £, with N^ 1 «C i <C 1. Since the macroscopic 
dynamics described by the GP equation is affected only 
by the structure of the wave function on length scales of 
O(l) and 0(1/N), these mesoscopic excitations have no 
influence on the evolution of the condensate. 
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